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Abstract. We analyse a system of self-gravitating identical bosons by means 
of a semirelativistic Hamiltonian comprising the relativistic kinetic energies of 
the involved particles and added (instantaneous) Newtonian gravitational pair 
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relativistic stability for such boson systems against gravitational collapse and 
to sharpen the predictions for their maximum stable mass. 
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1. Introduction 

In this paper we study the implications of two aspects of relativistic bound systems: 
the Coulomb (or gravitational) one-body coupling limit, and the effective coupling 
enhancement induced in a system of many identical particles interacting pairwise. 
These two effects lead to the conclusion that a system of N identical particles 
interacting by attractive 1/r pair potentials becomes unstable if N is very large. Our 
principal goal is to sharpen previous bounds on the critical mass of such a system. 

Relativistic quantum-mechanical theories imply an upper limit on the strength of 
the coupling of a single particle bound by an attractive Coulomb potential. Thus for 
a Hydrogen-like one-particle system with mass m, and units such that h = c = 1, the 
upper limits to the allowed coupling v in the potential — v/r are, respectively, v < 1 
for the Dirac equation and v < \ for the Klein-Gordon equation. Meanwhile, for 
the semirelativistic Salpeter equation [1-3] with Hamiltonian h = \J p 2 + m 2 — v/r, 
Herbst [4] showed that for v < 2/ir the spectrum of h in [0, m) is discrete and, 
moreover, he found an explicit lower bound. In summary 

2 

h = \/ p 2 + m? — v/r > myl — (ttv/2) 2 , v < — . (1) 

Under the Schrodingcr equation, with one or more particles, there is no such coupling 
restriction; thus the existence of such a coupling-limit is essentially a relativistic 
phenomenon. The Salpeter Hamiltonian has eigenvalues that lie between the 
corresponding Schrodingcr and Klein-Gordon energies. Thus, in addition to exhibiting 
the relativistic coupling limit, within the allowed couplings, the Salpeter energies are 
intermediate between those of Schrodinger and Klein-Gordon. For example, for the 
one-body problem with mass m = 1 and the potential V(r) = —v/r, the three theories 
have ground-state eigenvalues that depend on the coupling v as shown in Fig. 1. 
The Salpeter result was obtained by the use of a scale-optimized trial function with 
coordinate expression 4>{r) = ce~ T l a ; it is known analytically that the exact Salpeter 
curve is bounded below by the Klein-Gordon results for v < 1. A brief review of 
aspects of Salpeter semirelativistic theory may be found in Ref. [5]. If we consider, as 
we do in this paper, a system of N identical particles interacting pairwise via attractive 
potentials of the form —v/r^, then the necessary permutation symmetry of the wave 
function effectively enhances the pairwise coupling by a factor of the order of N. This 
effect, which we shall soon make clear, is most pronounced in the case of bosons. 

Particle identity in quantum mechanics is so strong that, in a system of identical 
particles, the particles lose their individuality; they cannot be separately tracked. 
This is often helpful for many-particle theory since, when it comes to permutation 
symmetry, at most two of the many possible Young Tableaux need be considered; 
moreover, many quantities are necessarily equal on the average. We shall now make 
clear the notion of effective many-body enhancement of the pair couplings which we 
alluded to above. We do this in the context of the problem that is the main concern 
of the paper. One of the advantages of the Salpeter semirelativistic theory is that it 
accommodates a straightforward formulation of the many-body problem. We consider 
therefore a semirelativistic system of TV self-gravitating identical bosons of mass m and 
momenta p i , i = 1, 2, . . . , N. This can be described — in a Newtonian approximation, 
justified to some extent by the assumption of a weak gravitational field — by the 
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Figure 1. Ground-state energies E (in dimcnsionless units) for the potential 
V(r) = —v/r according to the Schrodinger, Salpeter, and Klein-Gordon theories. 
The Salpeter curve is a variational upper bound. 
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H 



N 

E 

i=l 



+ m 2 



N 

l = i<j 



K 



, K > 0, 



(2) 



where, in the Newtonian pair potential, the gravitational interaction strength 
(determined by the gravitational constant G and the particle mass m) has been 
encoded in the coupling parameter k :— Gra 2 . The pair distance between the 
interacting particles i and j is given by = \xi — xj\. If we consider expectations 
with respect to a normalized boson function then we immediately find that there 
is a relation between H and a scaled two-particle Hamiltonian, namely, (H) = (H 2 ), 
where 

N 



Ho 



^vl+m 2 + ^v 2 +m 2 -(N-l) 



ri2 



(3) 



This expectation equality arises because the necessary boson permutation symmetry 
of the the exact TV-body ground state ^ implies that the expectations of the TV kinetic- 
energy terms in H are the same; and similarly for the |7V(7V— 1) pair-potential terms. 
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For convenience we have collected these into N/2 times a two-body Hamiltonian; we 
have included the overall N/2 factor and written this scaled two-body Hamiltonian as 
H 2 . Thus with respect to the exact wavefunction W we may write E = (H) = (H 2 ), 
where E is the corresponding exact energy. If we denote by E 2 the bottom of the 
spectrum of the two-boson problem with Hamiltonian H 2 , then, since boson symmetry 
in only two particles is in general a weaker constraint than symmetry in all N particles, 
it follows that E > E 2 . Indeed, the dependence of '5 on the variables {X3 . . . xjv} that 
are not present in H 2 cannot cause (H 2 ) to fall bellow the bottom of the spectrum of 
H 2 - Thus E 2 provides a lower energy bound to E. We shall sometimes express this as 
the operator inequality H > H 2 . A corresponding upper bound E g may be found with 
the aid say of a normalized Gaussian trial function $ g : thus E < E g = (Q g , HQ g ). 
These energy bounds allow us to compute bounds on the critical mass M c , the largest 
allowed mass for such a bound system. In order to make this point clear and to fix ideas 
we shall now compute an explicit energy lower bound and from this a lower estimate to 
M c . We first have to solve the two-body problem represented by H 2 . If we consider for 
this problem new coordinates R = xi + X2, and r = xi — X2, then the corresponding 
momenta are related by pi = p + P and \) 2 = p P. If we introduce a vector k 
which is orthogonal to p and P, then we may consider the following application of the 
triangle inequality: 

2{p 2 + m 2 )i = |2p + 2mk| 

= |p + P + mk + p - P + 77ik| 

< |pi + mk| + \p 2 + mk|. 

From this inequality and ([3]) we conclude the following inequalities 

(N- 1)k" 



H> H 2 > N 



V? 



+ vnr — 



Consequently from ([T]) we have 



E > Nm 



1 



f(N- \)kix 
V 4 



2r 



> Nm 



1 



Nun 



We have replaced N — 1 by N merely for analytical convenience. Thus we obtain the 
following A-boson lower energy bound 

Nktt 



Am , 
E > t(l 

TTK 



t := 



< 1. 



(4) 



It turns out that if we maximize the right-hand side of ^ with respect to N, that 
is to say, with respect to the parameter t, the critical value of t is t = l/\/2, so 
that the Herbst coupling inequality is satisfied at the optimal point. Since mass 
and energy are identified in our units, and k = Gm 2 , we arrive at the bound 
M c > (2/7r)/Gm ~ 0.63662/Gm. This detailed calculation shows how an energy 
bound leads to an estimate for the critical mass M c . The principal goal of the article 
is to refine such estimates. Thus we have here an explicit example of the phenomenon 
under discussion: if m is the mass of an alpha particle, say, then M cannot be larger 
than the mass of a modest mountain (we shall present an upper bound to M c shortly) . 
No such possibility arises from the corresponding non-relativistic theory. 



Semirelativistic stability of N -boson systems bound by l/r^ pair potentials 5 
2. Reduction: 'equivalent' two-body problems 

The question of the implications of the necessary permutation symmetry of the 
states for systems composed of many identical particles is almost as old as quantum 
mechanics. There are a number of historical threads. Before going into the technical 
details of our problem, we shall briefly mention two of these. The reasoning leading 
to the two-particle Hamiltonian H 2 suggests that the energy depends on a reduced 
density matrix p(xi, x 2 , x^, x 2 ) obtained by integrating 

*(X1,X 2 ,X 3 , . . . ,XAr)*(xi,X 2 ,X 3 , . . . ,Xjv) 

over all the variables x, with i > 2. A question raised by this is, what are the necessary 
features of p which characterize it as having come from an TV-boson function 'J? This 
is called the iV-representability problem and goes back at least to the early papers of 
P.-O. Lowdin [6] and A. J. Coleman [7]: a summary of early work in this direction can 
be found in the introductory chapters of Ref. [8] by A. J. Coleman and M. Rosina. 
Much of the early work was concerned with atomic and chemical systems. Density- 
matrix many-body theory and the A^-representability problem are still active areas of 
research [9-11]. We note that, as in the previous paragraph, one can derive energy 
lower bounds without attempting to solve the ./V-representability problem generally. 

Another story concerns nuclear type systems, where all the particles enter the 
motion on an equal footing, and considerations of centre-of-mass motion become 
important. In order to make this point more explicit, the example of the harmonic 
oscillator is helpful. We consider briefly the non-relativistic Hamiltonian given by 

N 2 N 
i=0 l=i<j 

The earliest treatment we know of for this problem is by Houston [12] in 1935; a 
solution expressed more specifically useful for our purposes was found in 1953 by 
Post [13]; the solubility of the iV-body harmonic-oscillator problem is periodically 
rediscovered, with justifiable fresh enthusiasm. In units with U = 1 the bottom of the 
spectrum is given exactly by the expression £ho — 3 (AT — l)y/Nv/(2m). The exact 
ground-state wave function is a Gaussian in the N — 1 orthogonal relative coordinates. 
If the same reasoning we used to derive the semirelativistic operator bound H > H 2 
above is now applied to ffno, the resulting lower energy bound obtained is exactly 
given by El = £ho / V%- If, instead, the 'reduction' (to a two-body problem) is effected 
with Jacobi relative coordinates, one obtains a lower bound for the harmonic oscillator 
equal to £ho itself. This type of reduction has its own history. In is only possible 
to indicate a few key events of this story in the present short article. In 1933, just 
after the discovery of the neutron, physicists began to look at few-nucleon problems. 
An approach emerged called the 'equivalent two-body method'. It was initiated by 
Wigner [14] and employed by many researchers [15-20] and eventually found its way 
into the pages of Rosenfeld's book Nuclear Forces [21] in 1948. The idea was always 
the same, to replace the TV-body problem by a tractable two-body problem. In many 
instances the result yielded an energy lower bound, but this was unknown to the 
workers at the time. The first rigorous results for such problems came in 1956 when 
Post [22] used Jacobi relative coordinates to show that indeed a lower bound could be 
constructed. In 1962 Post [23] applied this bound to the gravitational problem with 
pair potentials of the form —v/r^; together with a Gaussian trial function, the energy 
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was determined to 18%. Rigorous energy bounds with the aid of Jacobi coordinates, 
and a discussion of the 'equivalent two-body method' may be found in a paper by Hall 
and Post [24] in 1967. Some similar lower-bound results were later obtained by Levy- 
Leblond [25] and Stenschke [26]. An independent review and a certain sharpening of 
results by R. N. Hill may be found in Ref. [27]. The two streams of activity intersected 
in a paper by Hall [28] who used a similar argument to that of Coleman [7] to obtain 
a lower bound for a fermion system in the centre-of-mass frame in terms of a sum over 
N — 1 reduced two-particle energies; the bound was also optimized over a certain set 
of allowed non-orthogonal relative coordinates. Later this lower-bound theory (with 
non-orthogonal relative coordinates) was extended to the excited states [29] . A variety 
of alternative lower-bound models and approaches have been developed, for example 
by Carr [30], Manning [31], and Balbutsev [32]. The non-relativistic lower bound for 
the ground-state energy is rediscovered from time to time, for example by Membrado 
et al. [33] and Basdevant et al. [34]. 

3. Semirelativistic gravitating bosons 

We must now return to our main problem, the application of these ideas to a 
semirelativistic many-body system. The complication that the permutation symmetry 
(in its spatial aspect) is expressed in the individual-particle coordinates whereas the 
wave function is expressed in relative coordinates remains, and is adjoined by a new 
difficulty, namely the non-locality of the semirelativistic kinetic-energy operator. Here 
the TV-body harmonic-oscillator problem is now no longer solved exactly by the lower 
bound, but with a finite error is less than 0.15% [35]. 

The Hamiltonian H, among others, has been adopted to investigate spherically 
symmetric and nonrotating configurations of purely gravitationally interacting bosons 
forming compact objects known as "boson stars" [36-40]. This operator H is composed 
of the relativistically correct expression for the kinetic energy of all the involved 
bosons and static potentials n/rij which describe the gravitational forces between 
these particles. Therefore, it is clearly not possible in this model to take into account 
retardation effects. In addition, it goes without saying that this approach also omits 
general-relativistic effects [41-43]. Sufficient conditions have been found both for 
relativistic stability, which is characterized by the existence of a lower bound on the 
Hamiltonian H of ^ , and for relativistic gravitational collapse, which is inevitable 
if H is not bounded below. Moreover, semirelativistic bounds have been derived for 
the maximum possible, or critical, mass M c of boson stars, that is, the mass beyond 
which there must be relativistic collapse. 

The results of particular interest for this analysis can be summarized as follows. 
The relativistic kinetic energy \/p 2 + m 2 satisfies [36] a (tangential [44]) operator 
inequality, involving an arbitrary real parameter fi with the dimension of mass: 



This inequality can be adopted to relate the semirelativistic Hamiltonian H, Eq. ([2]), 
to its nonrelativistic counterpart. A variational bound on the ground-state energy 
of the nonrelativistic A^-particlc problem therefore translates into the upper bound 
M c < 1.52/Gm [37]. Exploiting the (only numerically computed) nonrelativistic 
ground-state energy, this bound is refined to M c < 1.518/Gm [39]. Rewriting H 
as a sum of one-particle Hamiltonians, each of which is bounded from below by the 



\/p 2 + m 2 < 



p 2 + to 2 + fi 2 
2/1 



V/i>0. 
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lowest positive eigenvalue of the Klein-Gordon Schrodinger equation with Coulomb 
potential, yields a bounaJH to the bottom E of the spectrum of H [38] : 




E ~V"'\/— y~ (N-l)n<l. (6) 

The replacement of TV — 1 under the square root by TV slightly weakens this bound 
but allows for its analytic maximization, which entails the analytic lower bound [38] 

4 0.7698 



M r > 



3V3Gm Gm 

Together these estimates constrain M c to the range 0.7698 < GmM c < 1.518. The 
resulting ratio of upper to lower bounds on M c is ry /l — 2.0. The so-called local-energy 
theorem may be used to increase the lower bound, whereas a more sophisticated choice 
of trial functions diminishes the variational upper bound. The combined effect of these 
improvements is to narrow down the range for M c to 0.8468 < GmM c < 1.439, with 
upper- to lower-bound ratio of ry/L — 1-7 [40]. 

We have re-analysed the upper bound of Ref. [40] with positive non-monotone 
Hartree wave- function factors <t>{r). With the factor (before scale optimization) <j)(r) = 
ce _r (l + ar), a > 0, we confirm the findings [40] that the best value of a is about 
a ~ 1, which yields GmM c < 1.43871 « 1.439. This <f> does indeed seem to be close 
to the best possible Hartree factor. With a = 1.13, we get a slight improvement, viz., 
GmM c < 1.43854. With the factor 0(r) = ce~ r (l — be~ r ) we obtain our best result, 
viz., GmM c < 1.43764 for b — 0.625. Thus we have been able to lower the upper 
bound on the critical mass slightly to M c < 1.438/Gm. 

In this paper we tighten the interval allowed for the critical mass of boson stars, 
by employing an improved analytic lower bound [45] on the ground-state energy of 
the TV-particle Hamiltonian ([2]) for semirelativistic self-gravitating TV-boson systems, 
to a range characterized by an upper- to lower-bound ratio ru/L — 1-3. The region 
of validity of a lower bound on the Hamiltonian H defines the range of relativistic 
stability of the gravitating TV-particle system under study [38]: our improved lower 
energy bound discussed below increases somewhat the stability region obtained in 
Ref. [38]; for instance, for couplings k < 1 — which allows for large TV — this increase 
of the stability range amounts to an 11% improvement. 



4. Lower bound for self-gravitating semirelativistic TV-boson systems 



Let | "J), 1 = 1, represent the normalized ground state of H, corresponding 
to its lowest eigenvalue E = (^f\H\^!). Now, the bosonic nature of the identical 
bound-state constituents forces the eigenstates of H (i.e., their wave functions) to 
be totally symmetric under any permutation of the individual-particle coordinates 
{xi,X2, . . . ,xjy}- The boson permutation symmetry of reduces the TV-body 
problem posed by the Hamiltonian H to a constrained two-particle problem [35]: 

TV(TV- 1) 



E = < * 



TV 



Pn 



7« 



7 



(7) 



By use of permutation symmetry, Eq. ((7|) may be cast into the equivalent form 



E = ( * 



(8) 



| Ref. [45] discusses this simple "N/2 lower energy bound" for arbitrary potentials. 
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After removal of the center-of-mass momentum from p! and p 2 , this apparent two- 
particle problem reduces to a one-body problem in the relative coordinate and 
momentum of the particles 1, 2 for which the Klein-Gordon equation with Coulomb 
interaction gives a lower bound: this eventually yields the bound 

The lower bound ((6]), however, is dramatically improved [45] by the use of Jacobi 
relative coordinates. The transformation from a given set {xj, i — 1,2, ... ,N} 



of coordinates to another set {p k 



1,2,. 



, N} may be defined by a matrix 
B = (Bki): p = Bx. The orthogonality B^ 1 — B T of B is not mandatory but 
may prove to be convenient. The momenta {iti} conjugate to the {Pj} are then also 
determined by tv = (B~ 1 ) T p = Bp. The transformation to Jacobi relative coordinates 
is represented by an orthogonal matrix with the first row given by 

' V» = l,2,. 



Bu = 



N 



whereas, for all 2 < k < N, in the fcth row only the first k entries are nonzero: 

1 



Bki — 



y/Kk - 1) 



Bkk — - 
B kl =0 



k- 1 

~k~' 
V i = k ■ 



Vi = l,2,. 



l,k + 2,. 



1. 



, N. 



Evidently, by the definition of B its first row generates the usual center-of-mass 
variable p x , while its second row introduces a pair distance p 2 = (xi — x 2 )/v / 2- 

Any boson state |$) is symmetric under permutations of all individual-particle 
coordinates. However, a non-Gaussian boson state is not necessarily symmetric in the 
Jacobi relative coordinates. Nevertheless, as has been shown in App. A of Ref. [45], 
each such |$) satisfies, for all i, k > 2, the 7V-representability identities 

<$| Pl . Pfc |$) =a jfc (s|pg|$>, 

($|7Ti-7r fc |$) =a< fc ($|ir||$>. (9) 
Now, for the sake of notational simplicity, let us introduce some abbreviations: 



- 1 




a "73T 


/ N 


N 




V N-l' 


N - 


2 




In -2 


N- 


1' 




/ N 



These parameters A, a, b, c are, of course, related by a 2 + b 2 = 2 and 1 + c 2 = 2A. In 
terms of Jacobi relative coordinates, the expectation value ([7]) then becomes 



E 



N 



<27Tl 



7T7V 



7K 



JV-ll 



We assume that the eigenstate \9) depends on {p 2 , p 3 , . . . , p^} but not on p ± . A 
lemma shown in Ref. [46] allows us to remove the center-of-mass momentum 7Ti from 
the kinetic term. Thus the A^-body ground-state energy E simplifies to 



E = ( * 



N 



7 AC 



\ap N - bp 



JV-ll 



(10) 
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Focusing to the (N — 1, N) subsystem we introduce new coordinates {R, r} and 
their conjugate momenta {P,p}, by performing the coordinate transformation 

R ^ _ O ( PN ^ { P\_0( 7T N 



r / V Pn-i J VP/ 2 V *n-i 
The expectation value (fTU)) suggests the most favourable choice of the matrix O: 

0= [ b ? ] =O l \ O l O = O z = 2. 



a —b 

Upon this change of variables, the ground-state energy E of the Hamiltonian H is 
given by the expectation value E = (^\Ti.\^) of the two-particle Hamiltonian 



1 K 



H = NJ{p + cP) 2 + m 

r 

It may be proved that Ti. is bounded from below by the Hamiltonian entering in the 
expectation value on the right-hand side of Eq. ([5]) (see App. B of Rcf. [45]). 

For the new momenta P and p, the identities ([9]) translate into the constraints 

($ |P 2 | $) = ($ |p 2 | $) and ($|P-p|$)=0. 

Consequently, we have to look for the bottom of the spectrum of the constrained 
problem posed by the operator Ti in a domain T> restricted by these conditions: 

V = {|^> e L 2 (5R 6 ) : {if |P 2 | tp) = |p 2 | ip) , (^|P • p|^> = 0} . 

This bottom £ of the spectrum of H, of course, provides a lower bound to E: 

£7 = <¥|«|*)> inf (tp\H\<p)=S. 

l¥>)ex> 

Let |^) G V, — 1, be the eigenstate of the Hamiltonian Ti. corresponding to 

this lowest eigenvalue £ . The eigenvalue equation of H satisfied by \tp) reads 

Ny/ (p + cP) 2 + m 2 \xjj) = (£ + ^ 

By squaring this relation and remembering the constraints that define T> we get 



jV 2 m 2 = 47 ( tp 



2 k£ 7K 2 
P ~ "27 ~ I^ 2 " 



^ . (11) 



Now, by assumption, is the lowest eigenstate of TL but not necessarily of the one- 
particle Kratzer-type [47] operator in Eq. (jlip . According to the variational principle, 
the (well-known) lowest eigenvalue of this Kratzer-type Hamiltonian provides a lower 
bound on the expectation value in Eq. (|11|) . Solving the implicit inequality for £ yields 
a lower bound to £, and thus to E [45]; this lower bound is nothing but the lowest 
positive eigenvalue of the corresponding Klein-Gordon Schrodinger equation [48] for 
gravitational interaction of appropriate strength: 




E > Nm\ ■ — , 7K 2 < 1. (12) 



Our improved lower bound (|12p on the ground-state energy of any self-gravitating N- 
boson system is of the same form as the bound © but with 7 = N(N — l)/2 replacing 
(JV - l) 2 , which is favourable since N(N - l)/2 < (N - l) 2 for N > 2. 
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5. Semirelativistic stability and critical mass of boson stars 

Let us now analyze the implications of the improved lower energy bound (|12p for both 
stability against gravitational collapse and maximum mass of boson stars. 

The existence of a lower bound on the spectrum of the Hamiltonian operator H 
guarantees the stability of the self-gravitating boson system against relativistic 
gravitational collapse. The region of validity of such kind of lower energy bound 
delimits the stability range of the bound state described by H. By construction, our 
bound (jT3J) holds for all N satisfying N(N — 1)k 2 < 2. This stability region is larger 
than the one, (N — 1)k < 4/-7T, found in Ref. [38]. For large values of N, allowed 
for sufficiently small couplings k, this gain amounts to tt/2\/2 = 1.11. In terms of 
Newton's constant G and the particle mass m, a sufficient condition for relativistic 
stability thus is that the particle number N fulfils the constraint 

Following Ref. [38], in order to allow for a discussion by elementary methods, we 
weaken Eq. (fT5|) by replacing the exact N dependence 7 = N(N — l)/2 by N 2 /2: 




E > Nm\ — — ^ '—, Nk<V2. 



Evidently, even this weakened lower bound is still above the lower bound ([5]) for 
all N > 2 + y/2 ~ 3.41; for large N, the weaker bound approaches the exact one. 
The (single) maximum of this weakened lower bound is situated at the critical point 
N = 4/3k, which is, fortunately, in the interior of the region of validity of our lower 
bound (fl"2|) on the Hamiltonian H as N < \f2j k. This maximum thus constitutes the 
(improved) lower bound on the critical mass M c of boson stars 

M > 4y/2m = 4V2 = 1.08866 
c ~ 3V3k 'SVSGm Gm 

This lower M c bound is larger by exactly a factor y2 than the result of Ref. [38]. 
Combining the lower bound (| 13[) with the Rayleigh-Ritz upper bound of Ref. [40] 
tightens the (Newtonian-limit) prediction for M c to 1.08866 < GmM c < 1.439, 
reducing thus the ratio between upper and lower bounds on M c to r v / L ~ 1.3. 

In summary, with the aid of an improved lower bound [45] (based on the relative 
coordinates of the bound-state constituents) on the bottom of the spectrum of 
the semirelativistic TV-boson Hamiltonian ([2]) with gravitational interaction we have 
succeeded in enlarging the range of semirelativistic stability of boson stars and in 
halving the theoretical uncertainty in the maximum mass of boson stars. 
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